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Abstract 

The area of a spherical region can be easily measured by considering 
which sampling points of a lattice are located inside or outside the re- 
gion. This point-counting technique is frequently used for measuring the 
Earth coverage of satellite constellations, employing a latitude-longitude 
lattice. This paper analyzes the numerical errors of such measurements, 
and shows that they could be greatly reduced if the Fibonacci lattice were 
used instead. The latter is a mathematical idealization of natural patterns 
with optimal packing, where the area represented by each point is almost 
identical. Using the Fibonacci lattice would reduce the root mean squared 
error by at least 40%. If, as is commonly the case, around a million lat- 
tice points are used, the maximum error would be an order of magnitude 
smaller. 

Keywords: Spherical grid; Golden ratio; Equal-angle grid; Non-standard 
grid; Fibonacci grid; Phyllotaxis 



1 Introduction 

The area of a region is easy to measure, without explicitly considering its bound- 
aries, by determining which points of a lattice are inside or outside the region. 
This point-counting method is commonly applied to estimating areas on a plane 
[5J [331 [2JJ E] • A related issue is how to approximate the region boundaries 
from this sampling [4]. 

Point-counting on the sphere is commonly used for estimating the Earth 
coverage of satellite constellations [36] [52], the fraction of the Earth's surface 
efficiently seen by one or more satellites. In the simplest case, each satellite 
covers a circular region (cap), so the constellation covers a complex set of iso- 
lated and/or overlapping caps |36j . Similarly, some maps designed for global 
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Figure 1: Latitude longitude lattice (top) and Fibonacci lattice (bot- 
tom), with 1014 and 1001 points, respectively. Orthographic projections, cen- 
tred at the pole (left), latitude 45° (middle) and equator (right). In the Fibo- 
nacci lattice, the points are much more evenly spaced, and the axial anisotropy 
is much smaller. 

earthquake forecasting depict earthquake-prone regions as complex sets of up 
to tens of thousands of caps [U [35J 121] • To assess these forecasts it is necessary 
to measure the fraction of the Earth's area covered by these regions [55] , 

Analytical solutions exist if the spherical region has a known, regular or 
polygonal shape [351 [SJ [53]. The area of a set of caps on the sphere has a 
complex analytical solution [36j , which unfortunately does not indicate whether 
any particular location on the surface of the sphere is covered. 

The numerical error of point counting should ideally decrease rapidly as 
the lattice density increases. Numerous works deal with errors on the plane 
[5J [33J [2S1 [23 [3] • Some particular cases using latitude-longitude lattices were 
analyzed elsewhere [36 . In automated counting, the computation time is di- 
rectly proportional to the number of lattice points. Satellite constellation cov- 
erage [53] or the areas marked on rapidly-updated earthquake forecasting maps 
[21] need continuous monitoring, implying a trade-off between accuracy and 
computational load. Thus it is important to find lattices able to measure areas 
as efficiently as possible. 

On the plane, the regular hexagonal lattice provides optimal sampling 12J. 
On the sphere, it is impossible to arrange regularly more than 20 points (the 
vertices of a dodecahedron), and the optimum configuration of a large number 
of points is problem-specific [BUI 121 EE [25J . For optimal point-counting, the 
area represented by every point should be almost the same. 

Traditionally, the latitude-longitude lattice is used for measuring Earth cov- 
erage 36, 53, 22J. However, it is very inhomogeneous (Fig. [I]), requiring non- 
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uniform weighting of the point contributions. Also, its number of points is 
restricted by geometrical constraints. 

The Fibonacci lattice is a particularly appealing alternative [J5] Ei El 1231 
E3liaESlEZlEHlCSlEai2HlESlESl- Being easy to construct, it can have 
any odd number of points [68], and these are evenly distributed (Fig. [JJ with 
each point representing almost the same area. For the numerical integration of 
continuous functions on a sphere, it has distinct advantages over other lattices 



This paper demonstrates that for measuring the areas of spherical caps, the 
Fibonacci lattice is much more efficient than its latitude-longitude counterpart. 
After describing both lattices (Sects. [2] and [3]), it is explained how to use them 
for measuring cap areas (Sect. [JJ and how the error of this measurement is 
assessed (Sect. [5]). The error results obtained from an extensive Monte Carlo 
simulation are described, and to some extent explained analytically (Sect. EJ). 
The conclusions are set out in the final section. 

2 Latitude— Longitude Lattice 

The latitude-longitude lattice is the set of points located at the intersections 
of a grid of meridians and parallels, separated by equal angles of latitude and 
longitude (Fig. [JJ). This is the "latitude-longitude grid" [68l[75] or "equal-angle 
grid" [25]. The points concentrate towards the poles, due to the converging 
meridians, resulting in high anisotropy. 

The number of points, P, depends on the angular spacing, 5, between grid 
lines. Since S = 180°/fc with k = 1, 2, . . ., 



That is the number of meridians (2k) times the number of parallels (k— 1), plus 
the two poles. Frequently, to evaluate satellite coverage, [55j 5 = 0.25°, so more 
than a million points are used. 

3 Fibonacci Lattice 

The points of the Fibonacci lattice are arranged along a tightly wound generative 
spiral, with each point fitted into the largest gap between the previous points 
(Fig. 0). This spiral is not apparent (Fig. [JJ) because consecutive points are 
spaced far apart from each other. The most apparent spirals join the points 
closest to each other, and form crisscrossing sets |68j . The points are evenly 
spaced in a very isotropous way. 

The next subsections describe how to construct the lattice used in this paper, 
and the history of the Fibonacci lattice in various research fields. 

3.1 Lattice Construction 

The Fibonacci lattice is named after the Fibonacci ratio. The Fibonacci se- 
quence was discovered in ancient India [62[ [40] and rediscovered in the mid- 
dle ages by Leonardo Pisano, better known by his nickname Fibonacci |61j . 



US] 



P = 2k(k - 1) + 2. 



(1) 
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Figure 2: Generative spirals of a Fibonacci lattice with 21 points. 

The angle turn between consecutive points along a spiral is based on the golden 
ratio ((/>): either the golden angle, 360°</>~ 2 ~ 137.5° (first spiral, red), or its 
complementary, 360°</>~ 1 ~ 222.5° (second spiral, blue). No point is placed at 
the poles. Orthographic projection, centred at longitude 0°, latitude 45°. 

Each term of the sequence, from the third onwards, is the sum of the previ- 
ous two: 0, 1, 1, 2, 3, 5, 8, 13, 21, ... Given two consecutive terms, Ft and 
Fi + i, a Fibonacci ratio is Fi/Fi+%. As first proved by Robert Simson in 1753 
|73j . this quotient, as i — > oo, quickly approaches the golden ratio, defined as 
cj, = 1 + (j,- 1 = (1 + y/E)/2 ~ 1.618. 

The Fibonacci lattice differs from other spiral lattices on the sphere [721 EH1 
l57l [60l [TT1 [6l [30] in that the longitudinal turn between consecutive points along 
the generative spiral is the golden angle, 360°(1 — (p^ 1 ) = 360° <\T 2 ~ 137.5°, 
or its complementary, 36O°0 _1 ~ 222.5°. Some lattice versions replace <f> by its 
rational approximant, a Fibonacci ratio. 

The golden angle optimizes the packing efficiency of elements in spiral lattices 
[58l [59] . This is because the golden ratio is the most irrational number [73] , so 
periodicities or near-periodicities in the spiral arrangement are avoided, and 
clumping of the lattice points never occurs [551 1151 HH (Ml HE] ■ 

The lattice version used here [55J is probably the most homogeneous one. It 
is generated with a Fermat spiral (also known as the cyclotron spiral) [TDJ, H5J 
IT?] , which embraces an equal area per equal angle turn, so the area between 
consecutive sampling points, measured along the spiral, is always the same 
|68j . Also, its first and last points are offset from the poles, leading to a more 
homogeneous polar arrangement [55J [55J than in other versions [55J 221 [521 [25J 
155] . When a Fibonacci ratio is used [55J [321 H21 EZH , the number of lattice points 
is F + l, where F > 1 is a term of the Fibonacci sequence. The lattice used here 
is instead based on the golden ratio and can have any odd number of points. 

To elaborate the lattice [55J, let N be any natural number. Let the integer 
i range from — N to +N. The number of points is 



P = 2JV + 1, 



(2) 



and the spherical coordinates, in radians, of the ith point are: 




(3) 



Ion.; = 2iri(j) 



-l 



(4) 
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This pseudocode provides the geographical coordinates in degrees: 
For i = -N, (-N + 1), . . . , (N - 1), TV, Do { 

lat, = arcsin (^if+i) x i 80 "/ 71 " 
Ion; = mod(i, (j>) x 360°/</> 
If Ion, < -180°, then Ion, = 360° + Ion, 
If Ion, > 180°, then Ion, = Ion, - 360° 
} End Do 

Here, arcsin returns a value in radians, while mod(a;, y) returns the remainder 
when x is divided by y, removing the extra turns of the generative spiral. For 
example, mod(6, <j)) = 6 — 3 X (f>. The last two lines keep the longitude range 
from -180° to +180°. 

Every point of this lattice is located at a different latitude, providing a more 
efficient sampling than the latitude-longitude lattice. The middle point, i = 0, 
is placed at the equator (lato = and lono = 0). Each of the other points 
(latijlon^) with i / 0, has a centrosymmetric one with (— latj, — lon^). The 
lattice as a whole is not centrosymmetric. 

The longitudinal turn between consecutive lattice points along the spiral 
(Eq. 0]) is the complementary of the golden angle ESI EE] ■ To use the golden 
angle instead, we can substitute Eq. H] by 

lonj = -2iri(j)~ 2 . (5) 

With Eq. HI the spiral progresses eastwards, while the minus sign of Eq. El 
indicates a westward progression. 

A remark not found elsewhere is that the lattice points are placed at the 
intersections between these Fermat spirals of opposite chirality, except at the 
poles (Fig. [5]). For drawing the spirals, i is made continuous in Eqs.EMHandEl 
and ranges from (— N — 1/2) to (N + 1/2). The 1/2 term accounts for the polar 

offset [211 EE]. 

3.2 Lattice History 

The Fibonacci lattice is a mathematical idealization of patterns of repeated 
plant elements, such as rose petals, pineapple scales, or sunflower seeds (Fig. [3]). 
The study of these arrangements is known as Phyllotaxis [SH GUI H > 43J . The 
Bravais brothers [TUj were the first to describe them using a spiral lattice on a 
cylinder. They argued that the most common angle turn between consecutive 
elements along this spiral in plants is the golden angle. The latter provides 
optimum packing [58, 59J, maximizing the exposure to light, rain and insects for 
pollination [IB]. Structures in cells and viruses also follow this pattern [50]. In 
some experiments, elements are spontaneously ordered on roughly hemispherical 
Fibonacci lattices, because the system tends to minimize the strain energy [IS] 
or to avoid periodic organization [18j . 

Unwrapping the cylindrical Fibonacci lattice produces a flat one [TUJ [551 EE] , 
frequently used for numerical integration [771 EH EH El] • 

Projecting the cylindrical Fibonacci lattice to the sphere generates the spher- 
ical version [55J EH] • This can be generalized to arbitrary surfaces of revolution 
[59l IT7] . The first graphs of spherical Fibonacci lattices used, as here, the 
golden ratio and a Fermat spiral [TSJ EE [T7] . A version based on the Fibonacci 
ratio is used in the modelling of complex molecules [TTJ ESI H2] • In this case 
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Figure 3: Spherical Fibonacci lattices in Nature an oblique view of two 
Mammillaria solisioides. The elements of many plants form Fibonacci lattices, 
for example, the areoles (spine-bearing nodes) of Mammillaria cacti [231 143) . 




Figure 4: The Starshine-3 satellite had 1500 mirrors arranged on a 

spherical Fibonacci lattice. Picture by Michael A. Savell and Gayle R. Fullerton, 
taken while the satellite was being inspected by John Vasquez. Reproduced by 
courtesy of the U.S. Naval Research Laboratory. 
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[551 W2\ , a "+" sign in the formula for the longitude should be substituted by 
"x" (D. Svergun, personal communication, 2009). Versions with the golden 
ratio serve to simulate plants realistically [23] and to design golf balls [75] . The 
latter method was used by Douglas C. Winficld (B. Moore and D. C. Winficld, 
personal communication, 2009) in the Starshine-3 satellite (Fig. |4j). 

The spherical Fibonacci lattice is a highly efficient sampling scheme for inte- 
grating continuous functions [5H HH1 H5J , as was observed in magnetic resonance 
imaging [2] . It is also advantageous for providing grid nodes in global meteoro- 
logical models [1H1 E3 E7J [55] EH [55] . 

4 Area Measurement 

The area measurement starts by finding which points of the lattice are inside 
the considered region. This is expressed by a Boolean function fi, such that 
fi = 1 if the ith point is inside, and ft = otherwise [36 . If the region is a 
cap with radius r, it suffices to measure d, the shortest distance (great-circle 
distance) between the sampling point and the cap centre (Fig. [5]): 

, f 1 if d<r, ,„s 
/i = \0if d>r. (6) 




Figure 5: Measurement of cap area. The cap (left), centred at p, has 
height h and great-circle radius r, and is a region of a sphere with radius R. 
Placed at random (right), its area can be estimated by considering the lattice 
points it covers. 

Each lattice point must be assigned a weight, Wi, proportional to the area it 
represents. Then the estimate, A, of the region area A, is measured considering 
the sphere area, As, and summing the contribution of all P points of the lattice: 

A^A = As ^ fm . (7) 

Ei=i 

The weights depend on the lattice type, as described below. 

4.1 Weights in the Latitude— Longitude Lattice 

The weights should be inversely proportional to the point density, which here 
increases towards the poles (Fig. [1]). The linear spacing between parallels is 
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constant. The length of a parallel is 2-kR cos(lat), where R is the sphere radius. 
In any parallel there is the same number of lattice points (2k), so their density 
is inversely proportional to [351 EH] : 

Wi = cos(lati) (8) 
4.2 Weights in the Fibonacci lattice 

Thanks to the even distribution of points in this lattice, the same weight can 
be assumed for all of them [56 , 68 , namely 

Wi = 1. (9) 

Each point represents the area corresponding to its Voronoi cell (Thiessen 
polygon) . This is the region of positions closer to the corresponding lattice point 
than to any other [51] HH] • Using the exact area of each Voronoi cell as weight 
for its lattice point [5] would improve the area measurement only slightly. The 
average cell area equals As /P. Here, regardless of P, only the areas of less 
than ten cells, located at the polar regions, differ by more than ~ 2% from this 
value [68] . As P increases, proportionally fewer cells depart significantly from 
the average area. Unlike the latitude-longitude lattice, the homogeneity of the 
Fibonacci lattice improves with the number of points. 



5 Error Assessment 

This section details how to assess the error involved in measuring the area of 
spherical caps placed at random on the sphere. A good way to measure the 
homogeneity of a spherical lattice is to compare the proportion of lattice points 
in spherical regions with the normalized areas of the regions [13 . For this task, 
it is natural to use spherical caps [601 1141 [S], which moreover appear in the 
applications mentioned in the introduction. 
The area of a spherical cap (see Fig. [5]) is: 

A c = 2nRh = 2ttR 2 (l - cos . (10) 

The normalized cap area is: 

A c 1 - cos(r/R) 
F -As~- 2 ' 

where As = AttR 2 is the sphere area. 

The absolute error of measuring a single cap is the absolute difference be- 
tween the estimated fraction and the actual one: 



E 



A c I A s - F 



(12) 



This depends on the lattice type, the number of points, and the size and location 
of the cap. If E = 0, the cap gets its fair share of weighted lattice points. If 
Ac = 0, E = 0. A plane cuts the sphere into two complementary caps. For any 
cap, E is the same as for the complementary cap with area A' c = 1 — Ac, so it 
suffices to consider caps not larger than a hemisphere (Ac = As/ 2). 
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The error is characterized here numerically using a Monte Carlo method. In 
each particular realization, a cap is randomly placed. Every point of the sphere 
has the same probability of becoming the cap centre, p, with coordinates: 



180° 



arcsin(2X — 1) 



(13) 
(14) 



7T 



lon B = 360 o X-180°. 



Here X is a random number, chosen with uniform probability in the range [0, 1], 
independently for each equation. The area of this jth cap is estimated, and its 
corresponding error is Ej . 

The process is repeated for a total of n independently located caps of the 
same size, providing a sample of n values of E. The root mean squared error is 



The supremum error in Eq. [12] for caps of any size, for a given lattice type 
and P, using Wi = 1, is the "spherical cap discrepancy" [5D1 [T3J [5] . ^ cann °t 
be determined exactly with a Monte Carlo simulation because it might result 
from a cap size not used, or a location not sampled. However, it is unfortunately 
difficult to compute explicitly [14 . The maximum E measured for the Fibonacci 
lattice is a lower bound to its spherical cap discrepancy. 

6 Results 

This section details the maximum errors and root mean squared errors measured 
with the Monte Carlo simulation detailed in the previous section. 

Thirteen lattice configurations of each type, from P ~ 10 2 to P ~ 10 6 , were 
analyzed. The chosen values of P increase in logarithmic steps, as accurately as 
possible (P ~ 10 2 , 10 7/3 , 10 8/3 , 10 3 , . . .). It is impossible to use identical P for 
both lattices because P is odd in the Fibonacci lattice but even in the latitude- 
longitude lattice. Moreover, there are considerably fewer possible values of P 
in the latter. For each configuration, 200 different cap sizes were used: from 
A c = 0.0025A S to A c = 0.5A S in steps of 0.0025A S . To obtain smooth results, 
n = 60, 000 was chosen for each cap size and lattice configuration. 

In total, 312 million caps were measured (2 lattice types x 13 values of P x 
200 cap sizes x 60,000 caps). After optimization, the calculation took 43 days 
of CPU time using 2.8 GHz, 64-bit processors. 

The maximum E measured, for caps of any size and location, is represented 
in Fig. [5] For the Fibonacci lattice, it is much lower and decays faster than 
for the latitude-longitude lattice, despite the non-uniform weighting of points 
of the latter. 

The rmse depends on the lattice type, number of points and cap size, as 
shown in Fig. [7J Because of the symmetry of the latitude-longitude lattice, any 
hemispherical cap covers one half of the points of the lattice, and (thanks to the 
point weighs) the estimation is perfect (E = and rmse= 0). This exception 
aside, the rmse tends to increase with the cap area, in a non-trivial way, which 
is more complex for the latitude-longitude lattice than for the Fibonacci lattice. 




(15) 
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Figure 6: Maximum error measured for randomly located spherical caps of 
any size. For area measurements on the Earth, about a million points are fre- 
quently used |22j , for which the maximum error would be an order of magnitude 
smaller in the Fibonacci lattice. 



Figure[5](top) shows the maximum values of rmse of each curve. They follow 
parallel power laws: 

rmse max ~ fcP~ 3 / 4 , (16) 

with k ~ 0.505 for the latitude-longitude lattice, and k ~ 0.362 for the Fibonacci 
lattice. Interpolating for the same P, the rmse max would be 0.505/0.362 ~ 40% 
larger in the latitude-longitude lattice. 

6.1 Analytical Approach to the Root Mean Squared Error 

The scaling shown in Fig. [8] can be explained using arguments from similar 
problems in the plane [371 1311 132] . The number of points of a regular square 
lattice enclosed by a sufficiently smooth, closed curve placed at random can be 
expressed as [52] : 

P in = AM 2 + D, (17) 
where A is the area enclosed by the curve, M is the inverse to the lattice spacing, 



and D is the discrepancy. The nominal spacing of a spherical lattice is <J 'As/P 
[55] , Assuming that the Fibonacci lattice is regular enough, 




Substituting this into Eq. [T7] 



(18) 



Pin « P-r + D. (19) 
As 
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Normalized cap area 
(fraction of sphere area) 




Normalized cap area 
(fraction of sphere area) 



Figure 7: Root mean squared error. For randomly placed caps which 
occupy an area fraction given by the abscissas, the curve indicates the root mean 
squared error of the area measurement. Each curve corresponds to a different 
number of lattice points, labelled at its right. For similar lattice densities, the 
Fibonacci lattice provides smaller and more homogeneous errors. 
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100 1000 10000 100000 1000000 

Lattice points, P 




100 1000 10000 100000 1000000 

Effective lattice points, P e 

Figure 8: Scaling of the maximum root mean squared error. Each dot 
corresponds to the maximum value of a curve in Fig. [7J Top: For the same 
number of points, the values would be about 40% smaller for the Fibonacci 
lattice. Bottom: The latitude-longitude lattice is inefficient because the effective 
(weighted) number of points is smaller than the real one. Considering this fact, 
the results collapse to a single function. The power-law decay and its exponent 
agree with analytical scaling arguments (see text). 



Dividing all the terms of this equation by P, 



\D\ rms(D) . , 

E^U- ; rmse«— ±-L. (20) 
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In the planar case, the root mean squared of D, rms(D), is proportional to 
y/~M [3"Tl |3~T1 I3"2] . Extrapolating this fact, we obtain the scaling observed in the 
Monte Carlo simulation: 

rmseoc^^^ocP-3/4. (21) 
P P y ' 

In the latitude-longitude lattice, the same scaling holds if we consider its 
smaller sampling efficiency. The latter may be measured with the denominator 
of Eq. This is the number of points of a homogeneous lattice that would do 
the same work, and is up to ~ 36% smaller than P in the range of P considered 
here. If Fig. [5] is plotted using these effective points in the abscissas, the same 
fit suffices for both lattice types (Fig. [5] bottom). 

If the sampling points were placed at random, the rmse would decrease more 
slowly, as oc p- 1 / 2 [7]. 



7 Conclusions 

This paper analyzes the errors involved in measuring the areas of spherical caps 
using lattices of sampling points: the latitude-longitude lattice (classically used 
for this task), and a Fibonacci lattice [68] . The latter has low anisotropy, is 
easy to construct, and is shown to result from the intersection of two genera- 
tive spirals (Fig. A a review of the literature (Sect. 13. 2[) reveals successful 
applications of this spherical lattice since the 1980s. 

If the Fibonacci lattice were used instead of its latitude-longitude coun- 
terpart, the area measurement would be more efficient, allowing a significant 
reduction of the computation time. For approximately the same number of lat- 
tice points, the maximum root mean squared error would be about 40% smaller 
(Fig. [SJ. The maximum errors would be also smaller, and would decay faster 
with the number of points (Fig. [5]). If about a million points were used, as is 
commonly the case |22) . the maximum error would be an order of magnitude 
smaller (Fig. [6]). 

It is also found that the maximum root mean squared error obeys a single 
scaling relation when the sampling efficiency is taken into account (Fig. [5]). 
This is partially explained using arguments from similar problems on the plane 

The area estimate depends also on the orientation of the sampling lattice, 
especially if the latter has high anisotropy. The difference may be assessed [28] 
by rotating the lattice [26) . Such an issue is not relevant in the case analyzed 
here because the caps were placed at random with uniform probability over the 
spherical surface. 

Here the Earth's shape has been approximated by a sphere [35], adding 
a slightly higher error than other shape models [38j \W[ [63]. Assessing this 
difference may be a topic of future research. 
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